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SUMMARY 


A computer code for the analysis of reinforced concrete structures has 
been developed using NASTRAN as a basis. Nonlinear iteration procedures were 
developed for obtaining solutions with a wide variety of loading sequences. A 
direct access file system was used to save results at each load step to restart 
within the solution module for further analysis. A multi-nested looping capa- 
bility was implemented to control the iterations and change the loads. -The 
basis for the analysis is a set of multi-layer plate elements which allow local 
definition of materials and cracking properties. 


INTRODUCTION 


Faced with an ever-growing need to ensure environmental safety, the 
nuclear industry must develop analytical tools to determine the physical integ- 
rity of reactor structures. These structures, typically layers of concrete 
interspersed with reinforcing materials, must withstand severe emergency loading 
conditions. Resulting high stresses due to these loads are expected to cause 
the formation and propagation of cracks within the concrete. Once this process 
begins, the structural characteristics change. Thus, a nonlinear solution al- 
gorithm is required with the capability of subjecting the structure to an arbi- 
trary sequence of loadings — thermal, pressure, gravitational and equivalent 
dynamics. Also, the algorithm designed must be efficient for the solution of 
large-order and complex structural models. 

Although the present NASTRAN (NASA Structural Analysis) computer program 
system contains a limited elasto-plastic analysis capability in Rigid Format 6 
(Piecewise Linear Analysis) , it was unacceptable for simulating the three- 
dimensional concrete cracking phenomena. The present NASTRAN stress-strain 
rules for nonlinear materials do not allow for a decrease in stress with in- 
creasing strain. This effect is certainly evidenced by the behavior of con- 
crete where the internal loads disappear due to cracking. Furthermore, the 
present analysis is limited to a stepwise application of a single load vector 
and requires large amounts of computer time. 
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This paper presents the implementation in Level 15.5 NASTRAN of a non- 
linear concrete analysis capability. The effort, performed by Universal 
Analytics, Inc. (UAI) and Ebasco, involved the addition of nonlinear elements, 
new functional modules, a DMAP (Direct Matrix Abstraction Program) alter pack- 
age, and modifications to the NASTRAN executive system. This project illus- 
trates the practicality of modifying NASTRAN for a specialized problem in 
contrast to the development of general capabilities as implemented in the 
standard NASA versions. 

The basic computational methods were chosen to provide efficient solutions 
to a wide variety of concrete cracking problems. Multi-layer quadrilateral and 
triangular plate elements with independent layer-by-layer material definition 
were implemented to calculate both linear and nonlinear stiffness matrix and 
load effects. 

A nonlinear iteration system was implemented to obtain solutions for 
large-order problems from efficient vector processing techniques. The system, 
controlled by DMAP, is controlled by a new module which performs nonlinear 
loading functions, convergence testing, and automated property updates (for 
design optimization). 

For user convenience in the practical solution of complex structural prob- 
lems, several modifications were implemented in the NASTRAN Executive system. 

A restart capability for re-analysis of a problem starting at any previous 
solution was performed with the addition of a random access storage file inde- 
pendent of the NASTRAN file system. New DMAP language modifications, currently 
in Level 16 NASTRAN, were installed to provide data block control for multi- 
nested DMAP looping. 

In the following discussion, the physical characteristics of the expected 
structures are described, followed by separate descriptions of the major items. 


PHYSICAL PROBLEM 


Many types of composite steel/concrete construction methods are used in 
the fabrication of high-strength structures. The hypothetical example shown 
in figure 1 illustrates the nature of a containment vessel used in nuclear 
power plant construction. Each layer represents a different property which 
may contain one or more types of materials. For example, a layer of closely 
spaced thin rods is imbedded in a concrete matrix and the stress-strain proper- 
ties represent the effects of both materials as a composite. The combination 
of bending and membrane forces in the structure contribute nonuniform stresses 
through the depth of the plate and therefore require a layer-by-layer analysis 
of the nonlinear effects. 

The nonlinear effects of concrete cracking are significantly different 
physically from the effects of plastic strains in metal structures. Classical 
elasto-plast ic analysis methods use one or more plastic "flow rules" which 
assume a continuous transition between elastic and plastic behavior with 
stresses always increasing with increasing strains (refs. 1 and 2). However, 
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the concrete stresses can be assumed linear until cracking occurs, at which 
point the stresses normal to the crack must be completely removed, leading to 
problems with discontinuities in the general solution. 

The theoretical limits for the principal strains in concrete are shown in 
figure 2, along with the actual approximations used. Also shown in the stress- 
strain diagram for a principal strain. The stresses, and therefore the result- 
ing internal forces, will be discontinuous and may lead to difficulties in the 
solution algorithm. This problem is discussed in a subsequent section. 


REINFORCED CONCRETE ELEMENTS 


New triangular and quadrilateral plate elements were developed to repre- 
sent the layered properties of the reinforced concrete. The generalized dis- 
placements and element geometry are identical to the existing NASTRAN TRIAl and 
QUADl plate elements. However, the nonlinear, multi-layered properties required 
a reformulation of the basic equations to account for the resulting coupling 
between the in-plane membrane and bending forces. A basic approach similar to 
the method used in reference 3 was adapted for use with NASTRAN- type element 
formulations. 

Each layer of the new triangular element (TRCA) illustrated in figure 3 
corresponds to a unique concrete or steel property which may be defined inde- 
pendently. The grid point locations may be off-set from the element reference 
coordinate system and the mean bending plane of the element does not require 
user calculation. 

In existing NASTRAN plate elements the membrane and bending properties are 
calculated independently. However, the coupling terms occurring in the layered 
elements require their simultaneous calculations. For these new elements, the 
existing NASTRAN bending routines have been expanded to include in -plane dis- 
placements and strains for each component basic triangle. Since the triangle 
and quadrilateral elements in NASTRAN are comprised of component triangles, the 
combination process also was expanded to include the coupled in-plane terms. 

During solution, the basic element calculations are used for several pur- 
poses. These include: 

1. Generation of a new stiffness matrix K to represent the current 
cracked condition. 

2. Generation of loads due to thermal expansion with forces across 
cracked layers removed. 

3. Generation of corrective forces P for each element representing 
the change in forces due to new cracks at a particular deformed 
shape . 

4* Calculation of stresses in the cracked layers for output processing. 
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5. If requested, the thicknesses of cracked layers may be auto- 
matically modified to eliminate the cracking strains. This 
option is used for design optimization. 

The finite elements which may include the standard NASTRAN elements are 
assembled into stiffness and mass matrices. Multipoint and single-point con- 
straints may be used to define boundaries. The static solution to applied 
loads, enforced displacements and temperatures is described in the next section. 


SOLUTION METHOD 


With finite elements, the general form of the static nonlinear equation of 
state is 

{$(u,T)} = {Pg} (1) 

where {i|>} is the generalized load vector resulting from element stresses 
u are the set of node point displacements 
T is the temperature field 

{P^} is the vector of externally applied loads at the node points 

The load function vector, $, may be obtained from the total energy of the 
system, U, and Lagrange’s equation, or: 




3U 

9u-r 


( 2 ) 


For the concrete cracking problem, the system may be reduced to a bounded 
linear problem. At any cracking condition, i, the generalized loads are: 

= [KiHui) - {Qi(T)} (3) 

where [K^] is the stiffness matrix and {Qi} is the thermal load vector. 

To be valid, the displacements {u^} in the above equation must be compat- 
ible with the cracking state used in [K^] . To find the displacements, the fol- 
lowing procedure is used. Equation (1) may be so structured as to produce the 
iteration equation below: 

[KaHui+i) = {Pg} - {$(ui)} + [Ka]{ui} (4) 

where the index i indicates each step of the iteration and each vector {u^+j^} 
is a new estimate of the displacement. 
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Two choices for the matrix [Kg] are available in the iteration process. 

In the Newton-Raphson method (ref. 4) the matrix [K^] is updated at each step 
using the new displacement vector, In the so-called perturbation 

method, a constant stiffness matrix [K^] may be used, thereby eliminating the 
need to decompose a new matrix at each time step. The perturbation method is 
used initially. If the measured convergence rate is too slow, the matrix is 
updated to speed up the convergence. 

If [Kg] is the matrix corresponding to the Initial cracking condition {u^} , 
the vector {$} from equation (3) may be substituted into equation (4) to produce 
the equation: 


tKaHui+i> = {Pg} + {Qi(T)} + [Ka - KiHu^} (5) 

Note that the stiffness matrix [Kg,] is constant, and no matrix decompositions 
are required except for the first step. 

For calculation purposes, equation (5) is rearranged to produce the actual 
iteration equation: 


[Kal{ui+i} = {Pg} + {Qa> + (6) 

where 

{6Pi} = [K^ - K.]{u.} + {Qi - Q^} (7) 

Note that the variable quantity on the right-hand side is the "corrective load" 
This is a function only of the elements which have changed in stiffness 
between states a and i. An exact solution is obtained when a displacement 
vector is found which satisfies equation (1), or if 

In the above method, the initial stiffness matrix [K^] may be obtained 
from any previous cracking condition or may be updated in the iteration pro- 
cedure. The starting vector (i = a) must be chosen such that no additional 
cracks will be formed with respect to state a. For this reason, the initial 
starting vector (u^^ at i = a) is always null. 

The iteration procedure is illustrated in figure 4 for a single degree of 
freedom problem with several layers of steel and concrete. The sawtooth curve 
represents the force-displacement function, $. Each peak represents the limit 
of each layer and the dashed line represents the applied load. Note that 
several solutions exist (where the dotted line intersects the solid lines) . 

If the structure was initially uncracked, the iterations would follow the 
dotted lines. However, once a layer is cracked, the unloading curve follows a 
different path, i.e., corresponding to the stiffness matrix of the cracked 
structure. New applied loads will also use this path. 

Both convergence and divergence of the system may be tested after each 
iteration. The loading error vector {6^} is defined as: 


A QQ 

■-to-/ 


{6i> 

and from equations (3) and (7): 


{$i> - {Pe> 


( 8 ) 


{6^} = [K^Hui} - {Qi> - {Pg} 

= - {6Pi> + [KgHuj^} - {Qg} - {Pg} (9) 

Substituting equation (6) for the second term on the right side, we obtain: 

{6j^} = {6Pi_i} - {6Pi> (10) 

The energy error of the system is obtained from the vector products of dis- 
placement and loads. As shown in figure 4 (i = 3) , the load error may diverge 
for one step, yet the process may still converge. In practice, however, for 
large-order problems, the curves $ become more continuous and the temporary 
type of divergence is rarely encountered. 

Many of the above matrix operations could be processed by existing NASTRAN 
DMAP as with the similar Level 16 differential stiffness procedure. However, 
because nonlinear element routines may be designed to directly calculate the 
corrective loads defined in equation (7), and because most of the vector opera- 
tions may be performed more efficiently in core, the basic iterations are per- 
formed in a new NASTRAN module. Described below is the basic organization of 
the NASTRAN implementation. 


NASTRAN IMPLEMENTATION 


The implementation of the concrete analysis (CA) algorithm in NASTRAN re- 
quired only isolated interfaces to the program due to NASTRAN ’s modular design. 
The interfaces of the new elements were straightforward due to the organization 
of the element table processing. Most of the new solution iteration procedure 
code was isolated to six new modules which used many of the existing NASTRAN 
matrix subroutines. 

The cracking analysis solution algorithm was implemented using the NASTRAN 
DMAP (Direct Matrix Abstraction) control language. The Level 16 DMAP compiler 
(ref. 5), available in Level 16 NASTRAN, was used to take advantage of its im- 
proved capabilities and constructs not available in the previous compiler. The 
most important of these features are the utility module SWITCH and the freedom 
to use Data Blocks before they are defined. SWITCH allows Data Block names to 
be interchanged so that two successive iterates may be SWITCHed and a branch 
made back to the iteration module. 

The initial framework for the cracking analysis DMAP used the Piecewise 
Linear DMAP of Rigid Format 6. DMAP ALTERS were inserted as required to 
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implement the new technique. The overall logical flow of the DMAP, paying par- 
ticular attention to the relationship of the new modules and loop controls, is 
shown in figure 5. The basic computational procedures for the cracking analysis 
are carried out by the six new modules that were implemented. A brief descrip- 
tion of each appears below. 

1. Functional Module CASMG (Strain Matrix Generator) 

Generates the cracking analysis strain matrix which will 
transform grid point idsplacements into generalized element 
strain of the element. This strain matrix is stored to be used 
within the iteration loop to compute: 

a. Thermal loads 

b. Element stress 

c. Corrective load vectors 

2. Functional Module CAUTIL (Cracking Analysis Utilities) 

Performs several utility functions for the cracking analysis. 

These include: 

a. Extracting data from the CA0F (defined in a later 
section) for previous material states 

b. Appending interpolated element temperatures to the EST 

c. Reorganizing property change specification Bulk Data 

d. Assembling partitioning vectors for CAITER (described below) 

3. Functional Module CASMA (Element Stiffness Assembler) 

Calculates the variable stiffness matrix for the cracking 
analysis. This matrix is assembled for each element on the 
initial pass, and only for updated elements on subsequent passes. 

4. Functional Module CAECPT (Cracking Analysis ECPT) 

Generates the updated ECPTCA data block with all element 
property updates for the current iteration. The updated ECPTCA 
is also saved on the CA0F. In addition, if DIAG 25 is set, a 
summary of all elements that cracked in the last stiffness 
iteration will be printed. 

5. Functional Module CAITER (Crack Analysis Iteration) 

Iterates to an equilibrium state solution for the nonlinear 
cracking analysis. 

6. Functional Module CA0UTP (Cracking Analysis Output Processor) 

Formats and prints the output data for a cracking analysis, 
including stresses, forces and property updates. 

Module CAITER is the primary iteration module. It controls the corrective 
load iteration, property updates, convergence testing, and setting of param- 
eters to control DMAP execution. The logical flow is shown in figure 5. 


In order to implement the desired solution paths of the cracking analysis, 
a random access ’Cracking Analysis Operating File’ (CA0F) was integrated into 
the solution methodology. The standard NASTRAN data base, G1N0, could not 
adequately handle the required procedures. The CA0F is a stylized version of 
the Automated Multi-Stage Substructuring Operation File (S0F) (ref. 6) found in 
standard Level 16 NASTRAN, This file allows the user maximum versatility during 
solution, the capability to use previous subcases as initial conditions for 
subsequent cases, and the cost effectiveness of internal restart. 

Data for each cracking element and each subcase is stored on the CA0F 
until deleted by the user. For elements, the stiffness and strain matrices are 
stored on the CA0F to allow rapid assembly of system matrices similar to Level 
16 NASTRAN. For each subcase, the current element summary (EST) and connection 
(ECPT) tables are stored. Also stored are the most recent load and displace- 
ment vectors which may represent a final or some intermediate stage in the 
iteration cycle. Current iteration parameters are also saved at each step. In 
the event that convergence has not been obtained, or a system failure has 
occurred, these parameters may be used to reinitiate the solution. This par- 
ticular powerful tool allows an ’internal’ restart to be performed either at 
some stage in the solution DMAP or within the GAITER module itself, outside the 
standard NASTRAN restart framework. 


CONCLUSIONS 


The nonlinear concrete cracking analysis described above has been imple- 
mented in a Level 15.5 version of NASTRAN. Although extensive testing has been 
performed using moderate size problems, the program is still being tested with 
large-order problems to provide experimental, correlations. With the basic 
cracking analysis capabilities and the additional user conveniences for multiple 
cases, the program is expected to become a powerful tool for the analysis of 
reinforced concrete structures. 

Many of the program techniques used in the concrete analysis modifications 
may be adapted to other types of NASTRAN analyses. A capability to perform 
automatic internal restarts inside a module, which could be implemented in 
other NASTRAN modules, has been proven. New methods of processing nonlinear 
elements have been implemented to provide additional efficiencies. New, highly 
efficient iteration techniques have been installed which combine core-held 
vector iterations in a module with out-of-core matrix iterations using DMAP. 
These operations could be adapted to perform more general types of nonlinear 
structure analysis such as elasto-plastic analysis combined with crack elements 
and nonlinear geometry behavior. 
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Figure 1.- Sample multi-layer reinforced concrete. 







Figure 3.- Triangular crack analysis plate element. 
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Figure 5.- Logical flow for concrete analysis in NASTRAN 
















Figure 5.- Logical flow for concrete analysis in NASTRAN. (Cont'd) 



























